Quantum simulation of the Dirac equation 
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The Dirac equation is a cornerstone in the history of physics merging successfuUy quantum 
mechanics with special relativity, providing a natural description of the electron spin and predicting 
the existence of anti-matter 0]. Furthermore, it is able to reproduce accurately the spectrum of 
the hydrogen atom and its realm, relativistic quantum mechanics, is considered as the natural 
transition to quantum field theory. However, the Dirac equation also predicts some peculiar effects 
such as Klein's paradox [3| and Zitterbewegung, an unexpected quivering motion of a free relativistic 
quantum particle first examined by Schrodinger Q]. These and other predictions would be difficult 
to observe in real particles, while constituting key fundamental examples to understand relativistic 
quantum effects. Recent years have seen an increased interest in simulations of relativistic quantum 
effects in different physical setups [^, 0, @, H, EIli ISi where parameter tunability allows 
accessibility to different physical regimes. Here, we perform a proof-of-principle quantum simulation 
of the one-dimensional Dirac equation using a single trapped ion 0], which is set to behave as a 
free relativistic quantum particle. We measure as a function of time the particle position and study 
Zitterbewegung for different initial superpositions of positive and negative energy spinor states, as 
well as the cross-over from relativistic to nonrelativistic dynamics. The high level of control of 
trapped-ion experimental parameters makes it possible to simulate elegant textbook examples of 
relativistic quantum physics. 

PACS numbers: 



The Dirac equation for a spin-1/2 particle with rest 
mass m is given by: 

dtp 

ih— = {c a ■ p + Pmc'^)tp . (1) 

Here c is the speed of light, p is the momentum oper- 
ator, a and /3 are the Dirac matrices, which are usually 
given in terms of the Pauli matrices [H, and the wave- 
functions ip are 4-component spinors. A general Dirac 
spinor tp can be decomposed into parts with positive and 
negative energies E = zty^p'^c^ + m^c^. Zitterbewegung 
is understood as an interference effect between the pos- 
itive and negative energy parts of the spinor and does 
not appear for spinors that consist consist entirely out 
of positive (or negative) energy. Furthermore, it is only 
present when these parts have significant overlap in po- 
sition and momentum space and is therefore not a sus- 
tained effect under most circumstances For a free 
electron, the Dirac equation predicts the Zitterbewegung 
to have an amplitude of the order of the Compton wave- 
length RzB ~ 10~^^ m and a frequency of w^b ~ 10^^ Hz 
and has so far been experimentally inaccessible. The real 
existence of Zitterbewegung, in relativistic quantum me- 
chanics and in quantum field theory, has been a recurrent 
subject of discussion in the last years [H, [l3]- 
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Quantum simulation aims at simulating a quantum 
system with a controllable laboratory system underlying 
the same mathematical model. In this way, it is possible 
to simulate quantum systems that can neither be efR- 
ciently simulated on a classical computer [isj nor easily 
accessed experimentally, while allowing parameter tun- 
ability over a wide range. The difficulties in observing 
real quantum relativistic effects have sparked significant 
interest in performing quantum simulations of their dy- 
namics. Examples include black holes in Bose-Einstein 
condensates Q , the Unruh effect in trapped ions @ , and 
Zitterbewegung for massive fermions in solid state physics 
01 . Graphene is studied widely in connection to the Dirac 
equation 



Trapped ions are particularly interesting for the pur- 
pose of quantum simulation [l^, [23, [2l|, as they allow 
for exceptional control of experimental parameters, and 
initialization and read-out can be achieved with high fi- 
delity. Recently, for example, a proof-of-principle simula- 
tion of a quantum magnet has been performed [2^ using 
trapped ions. The full three-dimensional Dirac equation 
Hamiltonian can be simulated with lasers coupling to the 
three vibrational eigenmodes and the internal states of a 
single trapped ion [8|. The setup can be significantly 
simplified when simulating the Dirac equation in 1+1 di- 
mensions, while the most stunning features of the Dirac 
equation remain, such as Zitterbewegung and the Klein 
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paradox. In the Dirac equation in 1+1 dimensions, 
dtp 

there is only one motional degree of freedom and the 
spinor is encoded in two internal levels, related to positive 
and negative energy states [H] . For the velocity of the free 
Dirac particle we find dx/dt = [x,Hr)]/ifi = cOx hi the 
Heisenberg picture. For a massless particle \ax^Iiu\ — 
and hence Ox is a constant of motion. For a massive 
particle [cr^;, Hj^ ^ 0, and the evolution of the particle is 
described by 

x{t) = x{0) + pc^H^^t + i e(e2»^o*/'* - 1), (3) 

with ^ = ^hc{(Tx — pcHjj^) Hjy^ . The first two terms 
represent an evolution linear in time, as expected for a 
free particle, while the last oscillating term may induce 
Zitterbewegung. 

For the simulation, a single ^"^Ca^ ion is trapped in 
a linear Paul trap [2^ with trapping frequencies ujax = 
27rxl.36 MHz axially and uJrad = 27rx3 MHz in the ra- 
dial directions. Doppler cooling, optical pumping, and 
resolved sideband cooling on the Si^2 ^ ^5/2 transition 
in a magnetic field of 4 G prepare the ion in the axial mo- 
tional ground state and in the internal state |S'i/2j™ = 
1/2). A narrow linewidth laser at 729 nm couples the 
states |1) := \Si/2:m = 1/2) and |0) := \Dr,/2,m = 3/2) 
which we identify as our spinor states. A bichromatic 
light field resonant with the upper and lower axial mo- 
tional sidebands of the |0) |1) transition with appro- 
priately set phases and frequency realizes the Hamilto- 
nian 

Hd = 2ijAn(TxP + hVla^. (4) 

Here, A = ^fijlmiOax is the size of the ground-state 
wave function with m the ion's mass, not to be confused 
with the mass m of the simulated particle, 77 = 0.06 is 
the Lamb-Dicke parameter, p = ^ is the momentum 
operator and and a are the usual raising and lowering 
operators for the motional state along the axial direction. 
The first term in Eq. [4] describes a state-dependent mo- 
tional excitation with coupling strength rjQ, correspond- 
ing to a displacement of the ion's wave packet in the har- 
monic trap. The second term is equivalent to an optical 
Stark shift and arises when the bichromatic light field is 
detuned from resonance by 251. Equation ([4]) reduces to 
the 1-f 1 dimensional Dirac Hamiltonian if we make the 
identifications c := 2r]flA and mc^ hfl. The momen- 
tum and position of the Dirac particle are then mapped 
onto the corresponding quadratures of the trapped ion 
harmonic oscillator. 

In order to study relativistic effects such as Zitterbe- 
wegung it is necessary to measure the expectation value 
{x{t)) of the harmonic oscillator. It has been noted theo- 
retically that such expectation values could be measured 



using very short probe times, \yithout reconstructing the 
full quantum state [1, [23, [1^. The unitary operator 
U = exp{—iriflp(Jxxt / A) , with x = {a^ + a) A and probe 
Rabi frequency flp, transforms the observable cr^ such 
that its slope for short probe times is proportional to (x) 
(see methods). To measure (x) we have to prepare the 
internal state in an eigenstate of Cj,, apply the unitary 
operator U for short times and record the changing exci- 
tation. Since the Dirac Hamiltonian generally entangles 
the motional and internal state of the ion, we first trace 
over the internal state (see methods). Next, we prepare 
an eigenstate of ay and apply the Hamiltonian generat- 
ing U for up to 14 fis, with 1-2 ^s steps, with the probe 
Rabi frequency set to ^lp = 2ttx 13 kHz. The slopes 
were obtained by a linear fit each based on IC* to 3 • 10^ 
measurements. 

We simulate the Dirac equation by applying Hjj for 
varying amounts of time and for different particle masses. 
In the experiment we set fl = 27r x 68 kHz, corre- 
sponding to a simulated speed of hght c = 0.052 A/fis. 
The measured expectation values {x{t)) are shown in 
Fig. [1] for a particle initially prepared in the spinor state 

tp{x;t = 0) = (\/27r2A)~^e~4A^ (J) by sideband cooling 
and application of a 7r/2-pulse. Zitterbewegung appears 
for particles with non-zero mass, obtained by varying f2 
in the range of < 17 < 27rx 13 kHz by changing the 
detuning of the bichromatic lasers. 

We investigate the particle dynamics in the cross-over 
from relativistic to nonrelativistic dynamics. The data 
in Fig. [1] were fitted with a model function of the form 
{x{t)) = at + RzB sinuJzBt- The fitted amplitude Rzb 
and frequency loz b of the Zitterbewegung are shown in the 
inset of Fig. [TJ From these data it can be seen that the 
frequency grows linearly with increasing mass ojzb ~ 217, 
whereas the amplitude decreases as the mass is increased. 
Since the mass of the particle is increased but the mo- 
mentum and the simulated speed of light remain con- 
stant, the data in Fig. [l]show the crossover from the far 
relativistic to nonrelativistic limits. Hence, the data con- 
firm there is vanishing Zitterbewegung in both limits, as 
theoretically expected. In the far relativistic case this is 
because luzb goes to zero, and in the nonrelativistic case 
because Rzb vanishes. 

The tools used for simulating the Dirac equation can 
also be applied to precisely set the initial state of the sim- 
ulated particle. The particle in Fig. [2^ was given an av- 
erage initial momentum {p(t = 0)) = h/A by a displace- 
ment operation with the Hamiltonian H = hrjflaxx/A. 
The initial state of this particle is built up out of a posi- 
tive energy component with positive velocity and a neg- 
ative energy component with negative velocity [2^. As 
long as these parts overlap, Zitterbewegung is observed 
which dies out as the parts separate. Further informa- 
tion is obtained by a complete reconstruction of the prob- 
ability distribution |'0(x)p shown in Fig. [2|3. It is also 
possible to initialize the spinor in a pure negative or pos- 
itive energy spinor (see methods). In Fig. [S^, the time 
evolution {x{t)) of a negative energy spinor with average 
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FIG. 1: Expectation values {x{t)) for particles with dif- 
ferent masses. The linear curve (■) represents a massless 
particle (fi — 0) moving with the speed of light given by 
c = 2ri^lA = 0.052 A/^s for all curves. The other curves 
are for particles with increasing mass moving down from 
the linear curve. Their Compton wavelengths are given by 
Ac := 2rjnA/n = 5.4A (▼), 2.5A (♦), 1.2A (•) and 0.6A 
(a), respectively. The solid curves represent numerical sim- 
ulations. The figure shows Zitterbewegung for the crossover 
from the relativistic 2riQ 3> to the nonrelativistic limit 
2riQ <C Q. The error bars are obtained from a linear fit as- 
suming quantum projection noise. The inset shows fitted Zit- 
terbewegung amplitude Rzb (B) and frequency lozb (•) ver- 
sus the parameter Q/riCl (which is proportional to the mass). 
Error bars la. 
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FIG. 2: Zitterbewegung for a state with non-zero average mo- 
mentum, (a) Initially, Zitterbewegung appears due to inter- 
ference of positive and negative energy parts of the state. As 
these parts separate, the oscillatory motion fades away. The 
solid curve represents a numerical simulation, (b) Measured 
(filled areas) and numerically calculated (solid lines) proba- 
bility distributions |?/'(x)p at the times t = 0, 75 and 150 /is 
(as indicated by the arrows in (a)). The probability distribu- 
tion corresponding to the state |1) is inverted for clarity. The 
vertical solid line represents {%) as plotted in (a). The two 
dashed lines are the expectation values for the positive and 
negative energy parts of the spinor. Error bars Icr. 



momentum (p) = 2.2fi,/A is shown. The corresponding 
reconstructed probability distributions are displayed in 
Fig. and it can be seen that there is indeed no Zitter- 
bewegung or splitting of the wavefunction. 

We have implemented a proof-of-principle quantum op- 
tical simulation of a tunable relativistic quantum me- 
chanical system. We have demonstrated that the simu- 
lated one-dimensional Dirac dynamics for a free particle 
shows Zitterbewegung and several of its counterintuitive 
quantum relativistic features. A natural route for the 
near future will be to move theoretically and experimen- 
tally towards the simulation of dynamics that are impos- 
sible (or difficult) to calculate in real systems, such as 
in quantum chemistry (27j or quantized Dirac fields in a 
quantum field theory context We consider this ex- 
periment to be an important first step that will pave the 
way towards more complex quantum simulations. Fur- 
thermore, the mapping between quantum optical systems 
and relativistic quantum mechanics may be followed by 
further analogies between the Dirac dynamics and the 
Jaynes-Cummings model [1, [2^ [l^ and photonic or 
sonic analogies [30|. 



Methods 
1 

Measurement of (x) and |V'(x)p. In ion trap ex- 
periments the only observable that can directly be mea- 
sured by fluorescence detection is cr^. Additional laser 
pulses can be used to map other observables onto (Jz- In 
the experiment we apply a state-dependent displacement 
operation U = exp(— z/cicra;/2) to the quantum state p 
followed by a measurement bf which is equivalent to 
measuring the observable 

^(fc) = U^azU = cos{kx)az + sm{kx)cry. (5) 

Here, k = 2riflpt/A is proportional to the interaction 
time t. If the ion's internal initial state is the eigenstate 
of az belonging to eigenvalue -1-1, {A{k)) = {cos{kx)) 
and for the eigenstate of a.y belonging to eigenvalue -1-1, 
{A\k)) = (sin(fcx)). A Fourier transformation of these 
measurements yields the probability density {S{x~x)) in 
position space. 

Fol: the position operator we have that ■^{A{k))\^_^ cx 
{xdy). Measuring {x) thus requires the preparation of an 
eigenstate of (jy which however cannot be done directly 
when the motional state is entangled with the internal 
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To distinguish between populations in the state |0) 
and |1), when reconstructing (as shown in 

Figs, m and [3] a short (200 ^s) fluorescence detection 
was applied to measure the internal state. Only cases 
where |0) was measured (leaving the motional state 
unchanged as no photons were scattered) were used 
for the subsequent analysis. To reconstruct |'(/'(a;)P 
belonging to |1), a 7r-pulsc prior to the short detection 
was used to interchange the internal state populations. 



Constructing a pure negative energy spinor. A 

general spinor is built up out of positive and negative 
energy components E± = ±-y/ c^p^ + m'^c^ such that 
ip = P^tp + P~ip. In momentum space, the projection 
operators are given by 
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FIG. 3: Time evolution for a negative energy eigenstate with 
Ac = 1.2A. Laser pulses create a negative energy spinor 
with average momentum (p) = 2.2ft/A. The corresponding 
initial momentum distribution |V'(p)P is shown in the inset. 
The filled curves are data, the solid lines represent a numer- 
ical calculation. The curve in (a) shows no Zitterbewegung. 
The solid curve represents a numerical simulation. In (b) the 
measured probability distributions |i/)(x)p for three different 
evolution times are shown (as indicated by the arrows in (a)). 
There is no splitting of the wavefunction and the evolution 
and spreading is as intuitively expected for a free particle. 



state. To solve this problem we first incoherently re- 
combine the internal state in This is done by first 
shelving the population initially in |1) to I-D5/2, m = 5/2) 
using a Rapid Adiabatic Passage transfer (RAP). A sec- 
ond RAP transfers the population in |0) to |1). A 
100 /is laser pulse at 854 nm transfers the population in 
|^5/2i ^ = 5/2) to IP3/2, m = 3/2) from which it sponta- 
neously decays to |1). The transfer efficiency is > 99 %, 
limited by the small branching ratio to the £'3/2 state. In 
the transfer steps a probability exists that the motional 
state of the ion is changed. This probability is however 
very small due to the small Lamb-Dickc parameter but 
could be eliminated completely by a separate measure- 
ment of the motional states of the spinor states |1) and 
|0), at the expense of a longer data acquisition time. 



(6) 



Here, I2 is the 2x2 identity matrix. The spinor state 
in Fig. [3] was engineered backwards by projecting out 
the negative energy part of a wavepacket with average 
momentum (p) = 2.2?i/A and renormalizing the spinor. 

The complete sequence for approximating the negative 
energy state is conveniently described in the basis of the 
eigenstates |±)j^ of ctj,. After ground state cooling we pre- 
pare the state |+)j,, then we displace this state to an aver- 
age momentum state (p) = 2.2?i/A by the displacement 
Hamiltonian H ~ hrjflayX / A. Next, a far detuned laser 
pulse rotates the internal state to G.84|-f-),y -t- z0.53|— )y. 
The displacement Hamiltonian H = —hrjQayx/ A shifts 
these parts in opposite directions to create the required 
asymmetry between the average momenta of the compo- 
nents. A final 7r/2-pulse creates the state shown in Fig.|31 
This state has > 99% overlap with the desired negative 
energy state. 
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